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Abstract 

Traditionally, it has been held that a central characteristic of stem cells is their ability to divide asym- 
metrically Recent advances in inducible genetic labeling provided ample evidence that symmetric stem 
cell divisions play an important role in adult mammalian homeostasis. It is well understood that the two 
types of cell divisions differ in terms of the stem cells' flexibility to expand when needed. On the con- 
trary, the implications of symmetric and asymmetric divisions for mutation accumulation are still poorly 
understood. In this paper we study a stochastic model of a renewing tissue, and address the optimization 
problem of tissue architecture in the context of mutant production. Specifically, we study the process of 
tumor suppressor gene inactivation which usually takes place as a sequence of two consecutive " hits" , and 
which is one of the most common patterns in carcinogenesis. We compare and contrast symmetric and 
asymmetric (and mixed) stem cell divisions, and focus on the rate at which double-hit mutants are gener- 
ated. It turns out that symmetrically-dividing cells generate such mutants at a rate which is significantly 
lower than that of asymmetrically-dividing cells. This result holds whether single-hit (intermediate) mu- 
tants are disadvantageous, neutral, or advantageous. It is also independent on whether the carcinogenic 
double-hit mutants are produced only among the stem cells or also among more differentiated cells. We 
argue that symmetric stem cell divisions in mammals could be an adaptation which helps delay the onset 
of cancers. We further investigate the question of the optimal fraction of stem cells in the tissue, and 
quantify the contribution of non-stem cells in mutant production. Our work provides a hypothesis to 
explain the observation that in mammalian cells, symmetric patterns of stem cell division seem to be 
very common. 

Author Summary 

Most cells divide symmetrically, by creating two copies of themselves. Stem cells have been considered an 
exception from this rule, because of their ability to divide in an asymmetric fashion, to create one stem 
and one differentiated cell. This asymmetric division mechanism became associated with the very concept 
of "sternness". In contrast to the traditional view, recent evidence suggests that asymmetric divisions, 
although common in some tissues, are less common in others. There, stem cell divide in a symmetric 
fashion, by either creating two copies of themselves, or two differentiated cells. What are the advantages 
and disadvantages of different types of self-renewal? One answer can come from examining evolutionary 
consequences of various stem cell division patterns. By using a mathematical model of cellular turnover in 
tissues, we studied the ability of cells to produce double-hit mutants (such mutants are often responsible 
for cancer initiation). It turns out that symmetrically-dividing cells generate such mutants at a rate which 
can be orders of magnitude lower than that of asymmetrically-dividing cells. This argument provides an 
evolutionary hypothesis of why stem cells in mammals often divide symmetrically. 

Introduction 

The ability of stem cells to divide asymmetrically to produce one stem and one non-stem daughter cell is 
often considered to be one of the defining characteristics of sternness. On the other hand, there is ample 
evidence suggesting that adult stem cell can and do divide symmetrically p][2]- 

Two basic models of stem cell divisions are discussed in the literature, see figure [TJ The asymmetric 
model suggests that the homeostatic control of the stem cell pool is maintained at the level of single cells, 
whereby each stem cell produces a copy of itself plus one differentiated cell 3 6 . From the engineering 
prospective, this model has the advantage of keeping the stem cell population level steady. An obvious 
disadvantage is its inability to replenish the stem cell pool in case of injury. This problem is naturally 
solved by the symmetric model, which maintains homeostatic control at the population level, rather 
than at the individual cell level. There, stem cells are capable of two types of symmetric divisions: a 
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proliferation division resulting in the creation of two stem cells, and a differentiation division resulting 
in the creation of two differentiated cells [TUTU]. Differentiation/proliferation decisions are though to be 
under control of numerous signals emanating from the surrounding tissue and the stem cells themselves 

mm. 



Asymmetric Symmetric 
Division Division 




Figure f : Symmetric and asymmetric stem cell divisions. In the asymmetric division model, a 
stem cell produces one differentiated cell and one stem cell. In the symmetric division model, a stem cell 
produces two differentiated cells or two stem cells. 

Uncovering division patterns of stem cells has been subject of intense research in the last fifteen 
years. Some of the first quantification of the division strategies in vitro comes from the work of [3D] 
who tracked mcthylation patterns in the dividing cells of the colon crypts. The analysis of the complex 
methylation patterns revealed that crypts contain multiple stem cells that go through "bottlenecks" 
during the life of the organism, which suggests that symmetric divisions are part of the picture. Another 
piece of evidence comes from experiments with chimeric mice to determine the dynamics of polyclonality 
of crypts. Initially polyclonal crypts eventually become monoclonal, which suggests that symmetric 
divisions must occur [3T][32]. By means of radiotherapy-induced mutations, the study of [33] suggests 
that a significant fraction of the somatic mutations in human colon stem cells are lost within one year. 

An important advance in quantification of symmetric vs antisymmetric divisions became possible 
with the invention of inducible genetic labeling [34]. This technique provides access to lineage-tracing 
measurements, from which the fate of labeled cells and their clones can be tracked over time. By means 
of the quantitative analysis of long-term lineage-tracing data [TOIES] ' ^ nas been shown that the rate 
of stem cell replacement is comparable to the cell division rate, implying that symmetric cell divisions 
contribute significantly to stem cell homeostasis [36f37| . The paper by [38] provides a review of the recent 
evidence of symmetric divisions in mammalian intestinal stem cells, spermatogenesis [35] and epithelial 
tissues such as hair follicles [39] . 

These new findings reveal that contrary to the previous thinking, adult tissue stem cells are often 
lost (e.g. by differentiation) and replaced in a stochastic manner. This notion challenges the traditional 
concept of the stem cell as an immortal, slow-cycling, asymmetrically dividing cell [34] . 

In paper [38], an important question is raised: Why should mechanisms of tissue maintenance so often 
lean toward symmetric self-renewal? One answer comes from recognizing the ability of symmetrically- 
dividing stem cells to respond to injury thus ensuring a robust mechanism of tissue homeostasis. It 
however could be argued that the symmetric divisions are "switched on" in response to a sudden stem 
cell loss, and the asymmetric division strategy is employed in the course of normal homeostasis. 

In the present paper, we explore an alternative hypothesis, which gives an additional "reason" for the 
tissue architecture favoring symmetric divisions. As a starting point, we note that in both symmetric and 
asymmetric division types, a dysregulation may lead to the loss of homeostatic control and an unchecked 
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growth of cells. A disruption in the control of proliferation/differentiation decisions can tip the balance 
and lead to abnormal stem cell expansion |40j . It has also been shown that disruption of asymmetric 
divisions can be responsible for cancerous growth of undifferentiated cells |4"TM45] . 

Here, we examine the symmetric and asymmetric divisions in the context of producing mutations. 
Many cancerous transformations start off by an inactivation of a tumor-suppressor gene [46| . This is 
the famous two-hit process discovered by Knudson [47l|48j and studied by many laboratories as well as 
theoretically. We ask the following question: from the point of view of two-hit mutant generation, what 
type of stem cell divisions is advantageous for the organism? What frequency of symmetric vs asymmetric 
divisions can maximally delay the stochastic generation of a dangerous mutant? To this end, we consider 
a continuous range of strategies with mixed type divisions and explore how the frequency of symmetric 
vs asymmetric divisions affects the generation of mutations. 

In this paper, we use both numerical simulations and analytical methods to study symmetric and 
asymmetric stem cell divisions in the context of mutation production. Other theoreticians have explored 
stem cell dynamics by means of deterministic stem cell modeling and stochastic numerical simulations 
[49H65] . A great review of many modeling approaches is provided in |66j . |67j studied the dynamics of 
mutation spread in development and showed that susceptibility to late-life cancers may be influenced by 
somatic mutations that occur during early development. |68| considered a model of stem cell dynamics, 
and calculated the rates of stochastic elimination (or washing out) of mutants. In this model, stem 
cells can proliferate symmetrically and differentiation is decoupled from proliferation. [51] considered the 
question of mutation generation by stem cells and found that mutations that increase the probability 
of asymmetric replication can lead to rapid expansion of mutant stem cells in the absence of a selective 
fitness advantage. 

In the present paper, we concentrate on the optimization problem of tissue architecture in the context 
of delaying double-hit mutant production, and focus specifically on symmetric and asymmetric stem cell 
divisions. We consider a stochastic model of double-hit mutant generation, and ask several questions 
related to evolutionary dynamics of mutations. What type of divisions is optimal? What cell types 
contribute the most to double-hit mutant generation? What is the optimal fraction of stem cells that 
delays carcinogenesis? 

Results 
Set-up 

We consider a two-compartment, agent-based model of stem cells and transit-amplifying (TA) cells. The 
stem cells are capable of both symmetric and asymmetric divisions (see figure [1]). The relative proportion 
of symmetric divisions can vary and is denoted by the symbol a (see Table [TJ, where a = 1 means that all 
divisions are symmetrical, and a = means that stem cells only divide asymmetrically. The symmetric 
divisions can be of two types, proliferation and differentiation. The type of symmetric division is defined 
by a regulatory mechanism which assures an approximately constant level of stem cells (sec Methods). 
The total population (which includes both stem cells, S, and TA cells, D) is denoted by N = S + D. An 
important parameter is A = N/D, which defines the proportion of stem cells with respect to TA cells: 
S/D = A — 1. 

We assume that the non-stem cells die, and that all cell types have a chance to divide. Each time 
a division happens, there is a probability, u\, that one of the daughter cells is a one- hit mutant with 
fitness r (while the fitness of all wild-type cells is given by 1). The fitness parameter defines the relative 
probability of the given cell-type to be chosen for division. In this paper we consider a range of fitness 
values, r, such that the one-hit mutants can be disadvantageous compared to wild-type cells, neutral, or 
even slightly advantageous. When a one-hit mutant divides, it has the probability ui to give rise to a 
two-hit mutant. Two-hit mutants are transformed cells which have a potential to give rise to a cancerous 
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Tabic 1: Model parameters 



Notation 


Description 


S,D 
N = S+D > 1 
A = N/D > 1 

< a < 1 
r < 1, r > 1 
Ui,«2 < 1 
-Ro-s-2 


Number of stem and non-stem cells 
Total population size 
Inverse relative number of non-stem cells 
Proportion of symmetric stem cell divisions 
Fitness of one-hit mutants 
Mutation rates leading to the acquisition of first and second hits 

Tunneling rate 



Notations used in the text and their brief description. 



tissue transformation. 

The generation of two- hit mutants is normally considered to be a rate-limiting step in cancer initiation. 
Once such a mutant is produced, it may break down homeostatic control and result in a wave of clonal 
expansion, followed by further transformations. It is this first step, the creation of a double-hit mutant, 
that we focus on in this paper. We investigate how the timing of such a mutant production depends on 
the tissue architecture, and specifically, on the symmetry of stem cell divisions. 

In order to gain analytical insights, a slightly simplified stochastic process was considered (see the 
Methods Section) which gave predictions that are in excellent agreement with the computational model. 



Tunneling rates 



While the detailed temporal dynamics of double-mutant production is given in the Methods Section, here 
we present the results for the so-called "tunneling rates" - the rates at which the stem cell system of 
a given size produces double- hit mutants (assuming that one-hit mutants drift at relatively low levels). 
Denoting the tunneling rate as i?o->2 (where the subscript suggests that the system transfers from all 
wild-type, "zero-hit", state to a system containing two-hit mutants), we have 
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where quantities y and ?/» satisfy the system 



= 

= r(l - u 2 )y 2 + A(l - y) - ry. 
The time to produce double-hit mutants is distributed exponentially with the mean 

To- 1 



(1) 



(2) 
(3) 



0^2 



Formula ((T|) describes the generation of double-hit mutants in the stem cells (the first term on the right) 
and in TA cells (the second term of the right). Several limiting cases are presented in Table [2] and 
illustrated in figure [2] 

Predictions of formula ([1]), as well as the more precise equation (jlip . have been compared with 
stochastic numerical simulations, and found to be in excellent agreement with them, see below. 
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Table 2: Important limiting cases for the tunneling rate (formula (flj) ). 



Regime 


Description 


Conditions 


i - y* 


1-2/ 


(IA) 


r < A, symm+asymm 


a > i*2, A - r > y^ui, r < A 


1 2u 2 
V CT (!-' r / A ) 


A-r 


(IB) 


r A, symm+asymm 


cr > y^tj, A — r| < y^ 






(IC) 

(2A) 
(2B) 
(2C) 


r > A, symm+asymm 

r < A, asymm 
r w A, asymm 
r > A, asymm 


a > |A - r\, |A - r\ > y^Tti, r > A 

cr -C U2, A — r > 012, r < A 
cr < y/uZ, A - r| < y^i 
a < |A - r|, |A - r| > y^tj, r > A 


Via- 1 ) 

1 <r(l-r/A) 
2u 2 
1 CT 

1 " 

2(r/A-l) 


1 — - 

r 

r«2 
A-r 

1 — - 

r 



The notations for the six different regimes refer to figure [2] 



(a) (b) 




io-' lir 5 o.ooi o.i io-' io- 5 o.ooi o.i 

Frequency of symm divisions, o Frequency of symm divisions, a 



Figure 2: The six different approximation regimes (Table [2]) for solutions of system ([2ll3]). 

Plotted is the quantity (a) 1 — y* and (b) y* as a function of the frequency of symmetric divisions, cr, for 
three different values of A (solid lines), together with the approximations given by the formulas in Table 
[2j Approximations (1A), (1-B), and (1C) are best demonstrated in panel (a), where the quantity 1 — y* 
is plotted. Approximations (2 A), (2B), and (2C) are best demonstrated in panel (b), where the quantity 
y* is plotted. The other parameters are Ui — it2 = 10 -7 , r = 1.1. 

Double-hit mutants are produced slower under symmetric compared to asym- 
metric divisions. 

An important question is how the fraction of symmetric divisions (a) affects the rate of double-mutant 
production. We can see that the production of double- mutants by non-stem cells does not depend on cr, 
the frequency of symmetric divisions. On the other hand, the production by stem cells is crucially affected 
by this parameter. Our formulas show clearly that the rate of tunneling grows as a decreases, and it is 
the highest when cr = 0, the case of purely asymmetric divisions. This means that in order to minimize 
the rate of double-hit mutant formation, one needs to maximize the share of symmetric divisions. In 
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figure [3] we plot the quantity 

Rate of double-hit mutant production by stem cells, under symmetric divisions ^o*^2 i ( <J = 1) 
Rate of double- hit mutant production by stem cells, under asymmetric divisions R^E'i' 7 — 0) ' 

for different percentages of stem cells. We can see that for realistic ranges of the mutation rates, the 
difference is at least 10-fold, and can be as high as 10 4 -fold, with the symmetrically dividing stem cells 
producing double- hit mutants slower than asymmetrically dividing cells. 




Figure 3: The reduction in the rate of double mutant production in stem cells with symmet- 
ric divisions compared to stem cells with asymmetric divisions only. Plotted is the quantity 
in formula (j4|) as a function of the mutation rate, iti. The percentage of the stem cells in the whole 
population (S/N) is marked next to the lines. The other parameters are 1*2 = Ui/2, r = 1. 

Figure U compares the analytical findings for the double-hit mutant production dynamics with the 
numerical simulations. We ran the stochastic numerical model (see Methods) for a fixed number of time- 
steps, and recorded whether or not a double-hit mutant has been generated. Repeated implementation of 
this procedure produced a numerical approximation of the probability of double-hit mutant generation, 
which is plotted (together with the standard deviations) as a function of er, the probability of symmetric 
divisions, for three different values of A, which measures the fraction of stem cells. Clearly, the probability 
of mutant generation in the course of a given time-interval is a decaying function of a. 

Another result that follows from our computations is the comparison of the double-mutant production 
in a hierarchical (stem cells plus TA cells) model compared with the conventional, homogeneous model 
that has been extensively studied |70H73| . It turns out the hierarchical model with purely asymmetric 
divisions always produces mutants faster than the homogeneous model. For the hierarchical model with 
purely symmetric divisions the result depends on the fitness of one-hit mutants. For disadvantageous 
one- hit mutants whose fitness satisfies r < 1, |1 — r| 3> ^/u2, the hierarchical model with purely symmetric 
divisions produces double-mutants faster, and for neutral and advantageous mutants, it produces double- 
hit mutants slower than the homogeneous model. In figure |4] we can see that for r = 1 (neutral one-hit 
mutants), hierarchical models with a sufficiently large values of a arc characterized by slower double- hit 
mutant generation compared to the homogeneous model (the horizontal line). 
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Figure 4: The probability of double-hit mutant generation as a function of a, the probability 
of symmetric stem cell divisions. The results of numerical simulations are presented as points 
connected with dotted lines (standard deviations are included) . Analytical results are given by solid lines 
(formula (fTTj) 1 ). The horizontal line represents the calculations for the homogeneous model. We ran 10 
batches of 1000 runs. The parameters are r = 1.0, ui = 0.00001, u 2 = 0.002, N = 500. 



Figure O shows additional results of simulations (together with our analytical calculations) , where 
for three different values of r (one-hit mutant fitness) the probability of double-hit mutant generation 
is plotted as a function of A. The values A — > 1 corresponds to a vanishingly low fraction of stem cells 
in the system, while A = 2 corresponds to 50% of all cells being stem cells. We show purely symmetric 
(a = 1) and purely asymmetric (a = 0) cases. For fixed mutation rates and populations sizes, the 
homogeneous model is characterized by only one parameter, r, which is the fitness of one-hit mutants. 
The probability of double-hit mutant generation strongly depends on whether these intermediate mutants 
are disadvantageous (r < 1), neutral (r ~ 1), or advantageous (r > 1). In contrast to the homogeneous 
model, the hierarchical model contains two additional parameters, A (the ratio of TA cells and the total 
population) and a (the probability of symmetric divisions). We can see that these two parameters affect 
the probability of double-hit mutant generation at least as strongly as the fitness r does. The influence 
of a is clear: the more the fraction of symmetric divisions, the slower double-hit mutants are produced. 
Next, we examine the role of the stem cell to TA cell ratio. 

The optimal fraction of stem cells. 

Let us consider an optimization problem for the tissue design, with the goal to delay the production of 
double-hit mutants. What is the optimal fraction of stem cells that the population should maintain? 
Analysis of the tunneling rates for a hierarchical model with purely symmetric divisions suggests that the 
optimal fraction of stem cells depends on the fitness of the one-hit mutants. If the one-hit mutants are 
disadvantageous (r < 1, [ 1 — r*] <C y/u2), then the tunneling rate grows with the parameter A. In other 
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Figure 5: The probability of double-mutant generation as a function of A, the ratio of TA 
cells to the total number of cells. As in figure 21 the results of numerical simulations are presented 
as points connected with dotted lines (standard deviations are included), and the analytical results are 
given by solid lines (formula The horizontal lines represent the calculations for the homogeneous 

model. We ran 10 batches of 1000 runs. Plotted is the probability of double-mutant generation as a 
function of A, for purely symmetric (a = 1) and purely asymmetric (a — 0) models, for three different 
values of r. The parameters are u\ = 10~ 5 ,U2 = 10~ 4 , N = 1000. 



words, in order to minimize the rate of double- mutant production, one would need to keep the stem cell 
pool as small as possible. 

For neutral and advantageous intermediate mutants, where the symmetric division model gives rise 
to the lowest double-mutant production rate compared to the homogeneous model and the hierarchical 
model with asymmetric divisions, this rate is minimized for a particular fraction of stem cells. This 
fraction is defined by the mutation rate ui in the neutral case, and by the fitness of the intermediate 
mutants in the case of weakly advantageous mutants. For neutral one-hit mutants (|1 — r\ <C y/vv), the 
optimal value of A is given by 

X opt = 1 + 2ul /3 , (5) 
and for weakly advantageous mutants with l<r<A, |r — 1| 3> i/u2, we have 

Aopt = (6) 

For example, for the biologically most relevant case of neutral one-hit mutants, the optimal fraction of 
stem cells is approximately 1% of the total population, assuming u-i = 10 -7 . 

These results are illustrated in figure [6] In this plot, we can see for r = 0.8 the probability of having a 
doubly mutated cell (after a given time-span) is an increasing function of A, as predicted. For the case of 
r = 1, the numerical simulation in figure |5] shows that X op t ~ 1.1 (compared with X opt = 1.093 predicted 
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by formula (0). For the case r = 1.2, formula ^ gives X opt ~ 1.5, which approximately coincides with 
the numerical optimum. In the case of advantageous mutants however the minima of A are very shallow. 




Figure 6: The probability of double-hit mutant generation in the symmetric division model. 

The case of symmetrically dividing stem cells, same as in figure [5] 



Do mutations in TA cells produce double-mutants? 

Let us compare the relative contributions to the double-mutant production rate coming from stem cells 
and TA cells, equation (fTJ): 

r&s = ^ (i - £) a - »), ^ = 4) (1 " y)- ^ 

The contribution from the TA cells grows as the fraction of TA cells increases. In figure [7] we plot the 
fraction of stem cells (given by 1 — 1/A) that corresponds to -Rq^.™ = #0^2 ■ We can see that for the 
mutation rates around 10 -7 , this fraction is about 0.1% for disadvantageous intermediate mutants, about 
0.5% for neutral mutants, and about 15% for advantageous mutants. This means that as long as the 
fraction of stem cells in the population is lower than these threshold values, TA cells contribute more to 
the production of double-hit mutants than stem cells. This threshold fraction grows for larger mutation 
rates, making it easier for TA cells to contribute significantly to the double-hit mutant production. An 
analytical approximation for the threshold value of A can be found for small values of mutation rates, 
such as 

A f 1 + ry^, r < 1, Regime (1A), 

\ r- ^f-, r > 1, Regime (1C). 

Next we address the question of optimization assuming that only mutations acquired by stem cells are 
dangerous and can lead to further malignant transformations. In this case, the rate of mutant production 



10 




Figure 7: The threshold fraction of stem cells corresponding to stem and TA cells con- 
tributing equally to double-hit mutant production. The quantity 1 — 1/A C , which corresponds to 
jjstem _ rTA^^ j g p^ted as a function of the mutation rate, U2, for three different values of r, and a = 1. 
For the fraction of stem cells above these values, stem cells have a higher contribution to the rate of 
double-mutant production compared to the non-stem cells. Thin dashed lines show the approximations 
of equation flS). 

is given by i?o^ ? 2 > equation ([?])■ It is easy to show that this quantity is maximized by asymmetric divisions 
only (a = 0), and it is minimized by symmetric divisions of stem cells (o = 1). Thus the message of this 
paper does not change if only stem cell mutations arc assumed to contribute to carcinogenesis. 

Discussion 

In this paper we found that symmetrically dividing stem cells arc characterized by a significantly lower rate 
of two-hit mutant generation, compared to asymmetrically-dividing cells. This is especially important 
in the context of tumor-suppressor gene inactivation, which is one of the more common patterns of 
carcinogenesis. This provides an evolutionary framework for reasoning about stem cell division patterns. 

In the literature, both types of stem cell divisions have been reported in various tissues. It has also 
been reported that the same stem cells are capable of both symmetric and asymmetric divisions. Whether 
a cell divides symmetrically or asymmetrically depends on factors such as the polarized organization of 
the dividing cell as well as the cell cycle length [73]. In Drosophila germ stem cells, cell division is 
asymmetric or symmetric depending on whether the orientation of the mitotic spindle is perpendicular 
or parallel to the interface between the stem cell and its niche [75) . Similarly, mammalian stem cells have 
been reported to employ both symmetric and asymmetric divisions to regulate their numbers and tissue 
homeostasis |76H77j . A switch from a symmetric mode of divisions to the asymmetric model has also 
been reported to take place in development (see [73[79] m the context of Drosophila). 

The fact that the rate of double-hit mutant production is the lowest for symmetrically dividing cells 
does not in itself explain or predict any aspects of the tissue architecture. It however provides an 
alternative hypothesis for the observation that in mammalian tissues, symmetric patterns of stem cell 
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division seem to be very common. The force of selection that comes from the cancer-delaying effect 
of such an architecture can be thought to have helped shape the observed division patterns. On the 
other hand, in more primitive organisms such as Drosophila, asymmetric stem cell divisions seem to 
dominate adult homeostasis (following the predominantly symmetric division patterns of development). 
Since cancer delay does not provide an important selection mechanism in the context of Drosophila, we 
can argue that this could help explain the observed differences. 

Symmetric divisions can have a cancer-delaying effect 

The mathematical result obtained here is that symmetrically dividing cells appear to delay double-hit 
mutant production compared to an equivalent system with asymmetrically dividing stem cells. What 
is the intuition behind this finding? Double-mutants are generated by means of mutations that happen 
in singly-mutated cells. To understand this process, let us focus on the dynamics of single mutants. In 
particular, we concentrate on singly-mutated stem cells, because the fates of single mutations in TA cells 
are identical in the two models. What happens to a singly-mutated stem cell under the different division 
patterns? 

If stem cells divide asymmetrically, then a mutation acquired in a stem cell will remain in the system 
indefinitely, because at every cell division, a new copy of the mutant stem cell will be generated. On the 
other hand, a mutant stem cell generated under the symmetric division model has a very different and 
much less certain fate. Each division of a mutant stem cell can result either in (1) elimination of the 
mutation from the stem cell compartment as a result of a differentiation, or (2) creation of an additional 
mutant stem cell as a result of a proliferation event. Superficially, it might look like the two processes 
might balance each other out. This intuition is however misleading. A lineage of mutant stem cells 
starting from a single mutant stem cell is much more likely to die out than to persist and expand. In fact, 
only 1/K of all such lineages will expand to size K. Half of the lineages will differentiate out after the 
very first division. Statistically there will be occasional, rare long-lived lineages, but the vast majority 
will leave the stem cell compartment after a small number of divisions. The production of those "lucky" 
long-lived mutants is not enough to counter-balance the great majority of the dead-end lineages that 
quickly exit the stem cell compartment. This is illustrated in figure [8j which plots the "weight" (the net 
size of a lineage over time, T) of a typical symmetrically dividing mutant stem cell, X sym , divided by 
the weight of a typical asymmetrically dividing mutant stem cell, X asym . The latter quantity is simply 
given by T, and the former quantity is a stochastic variable. We can see that the weight of symmetrically 
dividing mutant lineages is always lower than that of asymmetrically dividing lineages, which means that 
the former will have a lower probability to produce double-mutants offspring. We conclude that the 
uncertainty of fate of single mutant stem cells is the reason for the statistically longer time it takes for 
the symmetrically dividing stem cell model to produce a double-hit mutant. 

Interestingly, the above argument can be made in a similar manner for disadvantageous, neutral, or 
advantageous mutants. In any of those cases, an asymmetrically dividing mutant stem cell remains in the 
population indefinitely In the model with symmetric divisions, whenever a mutant stem cell is chosen 
for division, its probability to proliferate is similar to its probability to differentiate (in order to keep 
the homeostasis), and this the dynamics of each lineage is independent of its fitness (except that the 
frequency of updates is determined by the fitness of mutants; this is why the fitness parameter r factors 
out of equations J2]) and (|T2]l). 

We note that the effect of double- hit mutant production delay caused by symmetric divisions compared 
to asymmetric divisions is very significant. The difference in the tunneling rate which characterizes the 
time-scale of the process can be as high as 1, 000- fold for tissues with 10% of stem cells and the mutation 
rate of 10 per gene per cell division. 
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Figure 8: Why are symmetrically dividing stem cells produce mutants slower? The weight of a 
typical symmetrically dividing mutant stem cell lineage, X sym , relative to the weight of an asymmetrically 
dividing mutant stem cell lineage, X asym = T, is plotted as a function of the number of stem cell divisions, 
T. Here, S = 20, ./V = 1000, and 20 batches of 10, 000 simulations were performed to calculate the mean 
and the standard deviation. 

Can TA cells create double-hit mutants? 

The model studied in this paper tracks single- and double-hit mutant production in both stem and 
TA cells. It is interesting to compare which mechanism (through stem cell single mutants or TA cell 
single mutants) contributes more to the double-mutant production? It turns out that as long as the 
fraction of stem cells is smaller than a threshold (or equivalcntly, if the fraction of the TA cells is larger 
than a threshold), non-stem cells contribute equally or more to the production of double- mutants. This 
threshold fraction depends on (1) the mutation rate and (2) the fitness of intermediate, one-hit mutants. 
For example, if the intermediate mutants are neutral and the mutation rate is 10 -7 per gene per cell- 
division, then the threshold fraction of stem cells is about 0.5% of the total population. In other words, 
mutations originating in non-stem cells are significant if stem cells comprise less than 0.5% of the total 
population. This number is much higher if the intermediate mutants are advantageous, or if the mutation 
rate responsible for the second hit is higher. For U2 = 10~ 3 , non-stem cells are the driving force behind 
double-mutant production as long as stem cells comprise less than about 10% of the total population. 
This scenario is realistic in the presence of genetic instability, where inactivation of a tumor suppressor 
gene is likely to occur through a small-scale mutation of the first copy of the gene followed by a loss of 
heterozygocity event inactivating the second copy. The latter can happen a rate as high as 10 -2 [80] . 

The arguments presented above clarify some aspects of the long-standing debate about the origins 
of cancer, see also |81j . It is sometime argued that TA cells are unimportant for cancer initiation, for 
the following (quantitative) reason, unrelated to biological evidence. Intuitively, it seems that double-hit 
mutants cannot be created among TA cells, because all one-hit mutants in the TA compartment will 
be washed away before they have a chance to acquire the second hit. As John Cairns writes, "...there 
are 256 exponentially multiplying cells that divide twice a day and are being replenished continually 
by the divisions of a single stem cell, none of these 256 cells will ever be separated from the stem cell 
by more than eight divisions, and the replication errors made in those eight divisions are destined, of 
course, to be discarded", |82j . The computations in this paper demonstrate that under some realistic 
parameter regimes, double-hit mutants can be created in the TA compartment, and TA cells statistically 
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can contribute equally or more to double-hit mutant production compared to stem cells. The simple 
reason for this is as follows. Even though TA cells are short-lived, and getting a second mutation in a 
singly-mutated TA cell is unlikely, there are many more TA cells than stem cells. The low chance of 
double-mutant generation in a single TA cells can be outweighed by the fact that TA cells are a large 
majority, and single probabilities add up to create a significant effect. 

Cancer stem cell hypothesis 

The question discussed above is purely mathematical, and deals with the simple possibility to acquire 
two hits in the TA compartment. A related biological question is whether mutations occurring in the TA 
compartment can lead to further carcinogenic transformations, which brings us to the cancer stem cell 
hypothesis [8"51l8"4"] . While the concept of the cancer stem cell remains controversial [851181)] , here we do not 
intend to argue for or against this theory. Moreover, we refrain from making specific interpretations of this 
theory with regards to the exact origins of cancer. It has been argued that there is a distinction between 
the broader concept of the cancer stem cell on the one hand, and the narrower concept of normal stem 
cell becoming cancerous [53] . While the cancer stem cell hypothesis states that cancer is maintained by a 
small fraction of cells with stem-like properties, without making a specific assumption of how those cells 
are generated, the more narrow theory argues that mutations generated among non-stem cells cannot be 
cancer-initiating, because (at least, some) cancers originate via the creation of a cancer stem cell, which 
is a modified stem cell that retains some characteristics of " sternness" . 

In the light of this latter hypothesis, let us analyze the process of double-hit mutant production that 
occurs via mutations in stem cells only. Will our results change if only stem cell mutations can lead to 
carcinogenic transformation? To accommodate this assumption in our model, we must only use the first 
term in equation ([I]). It turns out that in this case, the message remains exactly the same: symmetrically 
dividing stem cell systems are characterized by a slower production of double-hit mutants compared 
to asymmetrically dividing stem cells. The universality of this result is explained above: the fate of 
mutations originating in the differentiated compartment is identical under the two models, and the only 
difference comes from the fates of mutant stem cells. 

Stochastic tunneling in the context of hierarchical tissue architecture 

Our theoretical results on the rate of double-hit mutant formation provide a generalization of a number 
of previous papers that studied the process of stochastic tunneling. The concept of stochastic tunneling 
was introduced by [70(171] when studying the first step in colon cancer initiation, the inactivation of the 
tumor suppressor gene APC. The concept has later been investigated by several groups in the context of 
cancer initiation, escape dynamics [72] , and more broadly as a means of crossing an evolutionary valley 
by an evolving species [73]. The basic Moran process in a homogeneous tissue has been used as the 
underlying mathematical model. A spatial generalization for the tunneling rate was calculated in [87], 
and a generalization to a specific model of renewing epithelial tissue was given in [81[I88]. The present 
paper expands the notion of stochastic tunneling to tissues consisting of stem and differentiated cells, 
whose fate can vary and is governed by relatively complex rules. Formula ([T]) includes the basic tunneling 
law of [701171] as a special case, and provides a way to predict the rate of mutant generation based on the 
stem cell fraction, the mutant fitness, and the probability of symmetric vs asymmetric divisions. 

Finally, we emphasize some of the important simplifications used in the present model. We consid- 
ered a two-compartment (stem/TA) model where all non-stem cells were treated as a single type. Our 
numerical explorations suggest that the addition of more compartments does not change the message of 
the paper, that is, in the presence of more cell types, symmetric divisions continue to minimize the rate of 
double-mutant production. Further, the effect of the stem cell niche was modeled in a very basic manner, 
by assuming the existence of a stem cell compartment and a relatively tight regulation of differentiation 
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vs proliferation decisions. Future directions include the addition of a more detailed description of spatial 
interactions, and the inclusion of other cellular processes such as de-differentiation. 

Methods 

Numerical simulations 

A stochastic numerical simulation was set up according to the following generalized Moran (constant 
total population) process. The population consists of four types of cells: stem cells (wild-type, i*, and 
one-hit mutants, j'*), and TA cells (wild-type, i, and one-hit mutants, j). We have i + i* + j + j* = N, 
where TV is a constant total population size. The dynamics proceed as a sequence of updates. At each 
update, one TA cell is randomly removed from the population, and replaced with an offspring of another 
cell, thus keeping the total population size constant. 

The process of division is modeled as follows. All cells (stem or TA cells) have a probability to divide. 
A cell is chosen for division based on its fitness. The fitness of mutated cells is given by r and the fitness 
of wild-type cells is 1. Let us use the notation J\f = i + i* + r(j + j*). Then the probability that a 
wild- type stem cell is chosen for division is given by i*/Af; the probability that a mutated stem cell is 
chosen for division is given by rj^/Af; the probability that a wild- type TA cell is chosen for division is 
given by i/Af; and the probability that a mutant TA cell is chosen for division is given by rj/Af. 

If a wild-type TA cell divides, it creates another wild- type TA cell with probability 1 — u%, and it 
creates a one-hit mutant TA cell with probability u\. If a mutant TA cell divides, it creates a one-hit 
mutant TA cell with probability 1 — U2, and it creates a two- hit mutant with probability ui. In case of 
such an event, the process stops. 

Divisions of stem cells can be either symmetric (with probability a) or asymmetric (with probability 
l — o"). Asymmetric divisions result in a creation of a TA cell. If a wild- type stem cell is dividing 
asymmetrically, then with probability 1 — u\ no mutations happen, and a one-hit mutant will be created 
with probability u%. In case of such an event, with probability 1/2 the TA daughter cell will get a 
mutation, and with probability 1/2 it will be the stem cell that acquires a mutation. Similarly, a one-hit 
mutant stem cell that divides symmetrically will create a two-hit mutant with probability it2, in which 
case the process stops. 

Symmetric divisions can be of two types: a differentiation, which results in a replacement of the 
dividing stem cell with two TA cells, or a proliferation which results in a creation of a stem cell. The 
probability of proliferation is taken to be p = gnS^^ J ! | _j- )io : where S is a constant parameter which 
measures the expected number of stem cells in the system. The probability of proliferation is given by 
1 — p. Again, when a wild-type stem cell divides, with probability 1 — u\ both daughter cells are wild- 
type, and with probability u\ one of the daughter cells is a one-hit mutant. If a one-hit mutant stem 
cell divides, both daughter cells are one- hit mutants with probability 1 — «2, and with probability u-^ the 
process stops because a double-hit mutant is created. 

The decision trees for stem cells are shown in figure HI for wild-type stem cells (a) and for mutated 
stem cells (b) . Stem cells are denoted by light circles with "S" and TA cells by shaded circles with "D" . 
One-hit mutants are marked with a star. 

These updates were performed repeatedly until either a double-hit mutant was created, or the maxi- 
mum number of time-steps was reached, which was set to 1000. We ran this code for 1000 times. After 
that we calculated the fraction of runs that resulted in a double- hit mutant, which approximates the 
probability of double- mutant creation. This quantity was calculated 10 times, and then the averages and 
standard deviations were calculated. 

To simulate the homogeneous Moran process, the same updates were performed except the number 
of stem cells was zero, + = 0. 
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Figure 9: Stem cell division decision trees for the numerical algorithm, (a) Divisions of wild- 
type stem cells, (b) Divisions of mutant stem cells. Stem cells are denoted by light circles with an "S" 
and TA cells by shaded circles with a "D" . One-hit mutants are marked with a star. 



Analytical tools 

Suppose we have the following version of the Moran process, which consists of a sequence of elementary 
updates. At each update, a daughter cell is chosen for death at random. Then a cell (a stem cell or 
a differentiated cell) is chosen to divide, according to its fitness, with mutants having fitness r. If a 
differentiated cell is chosen for division, it divides and this concludes the update. If however a stem 
cell is chosen for division, we proceed as follows. (1) With probability 1 — a, the stem cell can divide 
asymmetrically, which concludes this step. (2) With probability a, the stem cell divides symmetrically 
by differentiation, which is followed by a proliferation of another randomly chosen stem cell. Finally, 
another daughter cell is chosen for death, which concludes this step. 

The process described above is slightly different from the numerical agent-based algorithm outlined 
used in numerical simulations. In the generalized Moran process described here, the numbers of stem 
cells (S) and differentiated cells (D) are kept constant at every step. This is a simplification that allowed 
for analytical tractability (see below). In the numerical simulations the number of stem and differentiated 
cells fluctuates around a mean value, but despite this difference, the analytical formulas derived here are 
in an excellent agreement with the simulations. 

Note that in order to keep S constant, the symmetric stem cell divisions have to come in pairs (one 
proliferation and one differentiation event), and must be combined with two cell death events. Therefore, 
on the biological time-scale, an update involving symmetric divisions must have an average duration of 
two (and not one) elementary updates. Therefore below, when calculating various transition probabilities, 
the terms associated with symmetric divisions require a factor 1/2. 
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Let us denote by j* the number of single-mutant stem cells and by j the number of single-mutant 
differentiated cells. The updates can be envisaged as a Markov process in the space where 
J*, J > 0, with an additional state E denoting the generation of a double- mutant cell. Below we will 
use the condition that mutants are drifting at low numbers, j* S and j <C D. We have the following 
probabilities: 

• The probability that the number of mutant differentiated cells increases by one can be approximated 
as follows: 

r{j + (1 - a)j*) f D S ui miSct 

Pa a .a = (1 — U2l H + (1 — <7) , 

which is (i) the probability that a death of a wild- type differentiated cell (probability ~ 1), is followed 
by either a faithful division of a mutant differentiated cell, or a faithful asymmetric division of a 
mutant stem cell; (ii) a division of a wild-type differentiated cell with a mutation; (hi) an asymmetric 
division of a wild-type stem cell with a mutation happening in the differentiated daughter cell; (iv) 
a symmetric division of a wild-type stem cell with a mutation (times 1/2 by association with the 
symmetric division process). 

• The probability that the number of mutant differentiated cells decreases by one: 

j 



P 



D ' 



which is the probability that a mutant differentiated cell dies followed by a faithful division of any 
w.t. cell (« 1). 

• The probability that the number of mutant differentiated cells increases by two, and the number of 
mutant stem cells decreases by one: 

v - Gr i* ci ^ 

which is only possible for a symmetric update, when two w.t. differentiated cells die (probability 
« 1) followed by a mutant stem cell differentiating without a further mutation (probability ^-^(1 — 
M2)), followed by a w.t. stem cell proliferating without a mutation (probability w 1); the factor 1/2 
comes from the symmetric update. 

• The probability that the number of mutant stem cells increases by one: 

_ (1 - a)S ui a ( Su x rj* 

- n Y + 2 [~W + lv (1 ~ " 2 ' 

which is (i) the probability that following a death of a wild- type differentiated cell (sa 1), a wild- type 
stem cell divides asymmetrically with a mutation in the stem cell daughter cell, (ii) a wild-type 
stem cell proliferates with a mutation (ui), or (hi) a mutant stem cell proliferates without a further 
mutation (-^(1 — ^2))- 

• The probability to create a double-hit mutant: 

Pa a> e = \- (1 — a) h <y , 

which is (i) the probability that a mutant differentiated cell divides with a mutation, (ii) a mutant 
stem cell divides asymmetrically with a mutation, or (iii) a mutant stem cell undergoes either a 
differentiation or a proliferation event with a mutation. 
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Let us define by (fj r j (t) the probability to have j* mutated stem cells and j mutated differentiated cells 
at time t. The Kolmogorov forward equation for this function is given by 
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Let us define the probability generating function, 
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The probability to be in one of the states is given by ^(1, 1; t). Therefore, the probability to transit 

to state E is P2W = 1 — ^(L 1; t). The probability generating function satisfies the following first order 
PDE, derived by the standard methods (see e.g. [89]): 
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We have 



where 
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(l-y(f)) + |(l-y.(f)) 
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^-(y*+y 2 ) + K 1 -°')y*y (1-U2) 

K 1 - u z)y 2 + A(l - y) - ry, 

y(o) = 1. 



(10) 

(11) 

(12) 

(13) 
(14) 



Equation ([lip states that one-hit mutants in differentiated cells are produced by divisions of differentiated 
cells at the rate u\D and by divisions of stem cells at the rate u\S/2. The factor 1/2 comes from the fact 
that in asymmetric divisions, only a half of mutations will be in the differentiated cells, and in symmetric 
divisions which consist of pairs differentiation/proliferation, only half of the time a mutation will happen 
upon differentiation. Mutations in stem cells are produced by the divisions of stem cells at rate u\S/2. 
The ordinary differential equations describe the dynamics of lineages that start from one differentiated 
mutant (equation for y) or from one stem cell mutant (equation for y„). The dynamics of differentiated 
mutants is independent of o. 

Let us first solve equation (|13[) , which informs us about the probability of creating a double- hit mutant 
in a differentiated cell. This Riccati equation can be solved by standard methods, and the growth of the 
quantity 1 — y proceeds in the following stages: 



The linear growth stage, where 1 — y ps ru2<, as long as t <C f* (to be defined). 
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• The saturation stage, where 1 — y ps C , as long as t ^> i„. 
The constant C obtained from the stable fixed point of equation (|13[) is given by the equation 



_ A + r - y/(A + r) 2 - 4rA(l - u 2 ) 
V ~ 2r(l - u 2 ) 

and can be approximated by concise expressions as shown below. Given the solution for y, equation (JTJ 
can also be analyzed. The function 1 — y t increases monotonically and reaches saturation at 1 — y* = 1, 
after characteristic time £**. To find that time-scale, we substitute the constant approximation for the 
function y, to obtain 

y*(t) ps exp(—t/U*). 

There are several regimes where the expression take a particularly simple form (sec Tabic [TJ. 

Regime (2A). Let us assume that |A — r\ \fu2~, r < A, and a <C U2- In this case, we have 

ru2 1 A — r 

A — r A — r \r1i2 

There are therefore three distinct regimes defined by the behavior of the functions y(t) and y*(t). 

1. If t <C i*, we have 1 — y = r«2i and 1 — y* = ruit + (r£) 2 M2/2. In this case we have 

Nuiiizrt 2 Su\U2r 2 t i " 
2 12 - 



Pr(*) = l-esp - 



where the second term in the exponent is typically smaller than the first, and the behavior is thus 
indistinguishable for the usual homogeneous Moran process at early times. 

2. If t* -C t -c i**, we have 1 — y = ru2/(\ — r) and 1 — y* = \ru2tl (A — r). In this case we have 

P^-(i) = 1 - exp (J D + S / 2 )W* _ gWr^A . (15) 
\ A — r 4(A — r) J 

3. Finally, if t ^> we have 1 — y = ru2/(X — r) and 1 — = 1, and 

P-(t) = l-e^(-PW), Po^= (J + f /2)W + ^. (16) 

A — T z 

This regime becomes unimportant if for t = i** we can show that the quantity in the exponent is 
much larger than one. We have 

PT'(^) = Pf-(U) = 1 - exp f- (jP + S/2)U1 -< SUliX r) 



r \ru2 

and this quantity is very close to 1 for example if u\ ~ U2 and (A — r)S 1. 
Regime (2B). Let us assume that |A — r| -C \fu2~ and cr <C \fuz- In this case, we have 

C = yf^2~, U = ——, t** = -— 

2A^/u^ A^/m 2 

There are therefore only two regimes defined by the behavior of the functions y(t) and y*{t). 
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1. If t <C t*, we have as in the previous case, 1 — y = rw2t and 1 — y* = rit2i + (rt) 2 U2/2. The 
probability of double-hit mutant production is thus given by 

P?"(t) = 1 - exp _ **^!) i 

where the second term in the exponent is typically smaller than the first, and the behavior is thus 
indistinguishable for the usual homogeneous Moran process at early times. 

2. If <C t, we have 1 — y = y/u% and 1 — = 1. In this case we have 

P 2 sat {t) = l-exp{-R ^ 2 t) 1 R^ 2 = {D + S/2) Ul ^T 2 + ^. (17) 

Regime (1A). Let us assume that |A 
as a linear function, 



for t+ <C t < i**, where 
For t > the quantity 1 — y*(t) tends 

l -^ = {^hy < 19 > 

Note that the initial behavior of the function 1 — y* (t) does not depend on a. This means that for relatively 
short times (t <C i**), the mutant generation in stem cells proceeds in the same way for symmetric and 
asymmetric divisions. The length of this regime and the level of saturation however are both functions 
of a. It is easy to see that both t** and the saturation level increase as a decreases. This means that the 
rate of mutant accumulation becomes higher for asymmetric divisions. 

Regime (IB). Let us assume that |A — r\ <C \fu2 and o S> ^/u2- Now, the linear stage for 1 — y*{t) 
is defined as 

and it occurs for the times f,<fCf», where 

1 _ 1 

For < > f„. the quantity 1 — y*{t) tends to a constant, 

i-y.(*) = yf«a /4 . (20) 

Calculations for regimes (1C) and (2C) are performed in a similar manner, see Table [2J 



r\ yfu-i; r < A, and a 3> u%. The quantity 1 — y*(t) behaves 
, , \ru 2 t 

i-y*{t) = ^- J -, 

A — r 



1 1 / A — r 

t** = -\ r— . (18) 



A — r r V 2ii2Afi 

to a constant 
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